Setup

Load necessary packages and clear workspace

library(ggplot2)
library(viridis)
library(lme4)
library(lmerTest)
library(emmeans)
library(knitr)
library(dplyr)
library(drc)
library(aomisc)
library(nlme)
library("ggpubr")
library(broom)
library(gridExtra)
library(gtable)
library(grid)
library(psycho)
rm(list = ls())

Setup plot themes

#Set text size and colors
txt.size = 12
txt.color = "black"
txt.face = "bold"

###define specs for agreement values plots over sample sizes
#agreement value colors
d.color = "#FFC700" #actual data 
random.color = "#00B885" #random data
crossMov.color = "#4D4D4D" #agreement values when actual data compared to cross movie.
factor.order <- c('crossMov', 'random', 'sample')
factor.color <- c(crossMov.color, random.color, d.color)
factor.shape <- c(24, 23, 21) #triangle, diamond, circle
factor.label <- c('Cross-movie', 'Random', 'Same')
#jitter for individual agreement values
jitter.alpha = .1 
jitter.width = 0.2
#point for mean agreement value 
point.size = 5
point.alpha = 1
#errorbar width 
errorbar.width = .8
errorbar.alpha = 1
#elbow related plot 
elbow.size = 8
elbowPos.scale = .01

#define specs for curvature function fittings. 
vir_7 = viridis(n = 7) #load the first 8 viridis colors
line.colors <- c(vir_7[1],vir_7[1], vir_7[2],vir_7[3], vir_7[4], vir_7[5], vir_7[6]) #Linear, Exponential, logarithmic, Asymptotic, Power Curve, Yield loss.
dot.color <- vir_7[7] #color for dotplot 
fit.jitter.alpha = .05
fit.jitter.width = 0.3
fit.line.size = 1.5
fit.line.alpha = 1

###define specs for saving agreement and function fit plots
figure.width = 25
figure.height = 19 
figure.unit = "cm"
figure.bg = "transparent"
figure.filetype = "png"

grain.labs <- c("Coarse", "Fine")
names(grain.labs) <- c("c", "f")

sampleSize = seq(from = 2, to = 32, by = 2)

#Set theme for plotting 
theme.esMethods <- theme(
panel.grid.minor = element_blank(), 
panel.grid.major = element_blank(), 
#panel.background= element_blank(), 
panel.border = element_blank(),
panel.spacing = unit(.05,'in'), 
panel.background = element_rect(fill = "transparent",colour = NA),
plot.background = element_rect(fill = "transparent",colour = NA),
axis.line = element_line(size = .2, colour = "black"), 
axis.title = element_text(size = txt.size),
axis.text = element_text(size = txt.size),
strip.background = element_rect(colour = NA, fill = "transparent"), 
strip.text = element_text(size = txt.size), 
legend.title = element_blank(),
legend.key = element_rect(colour = NA, fill = NA),
plot.title = element_text(hjust = 0.5)
)

Define necessary functions

Function: elbow_finder Adapted from https://stackoverflow.com/questions/2018178/finding-the-best-trade-off-point-on-a-curve/2022348#2022348

Function description: to find the stabilization point (elbow) of change in segmentation agreement based on best fitting growth/ decay functions.

Returns: coordinates of curve that is furthest from the straight line connecting the minimum and maximum points.

Parameter x_values: array of x coordinates of curve

Precondition x_values: numeric

Parameter x_values: array of y coordnates of curve

Precondition y_values: numeric

elbow_finder <- function(x_values, y_values) {
  # Max values to create line
  max_x_x <- max(x_values)
  max_x_y <- y_values[which.max(x_values)]
  min_x_x <- min(x_values)
  min_x_y <- y_values[which.min(x_values)]
  max_df <- data.frame(x = c(min_x_x, max_x_x), y = c(min_x_y, max_x_y))

  # Creating straight line between the max values
  fit <- lm(max_df$y ~ max_df$x)

  # Distance from point to line
  distances <- c()
  for(i in 1:length(x_values)) {
    distances <- c(distances, abs(coef(fit)[2]*x_values[i] - y_values[i] + coef(fit)[1]) / sqrt(coef(fit)[2]^2 + 1^2))
  }

  # Max distance point
  x_max_dist <- x_values[which.max(distances)]
  y_max_dist <- y_values[which.max(distances)]

  return(c(x_max_dist, y_max_dist))
}

Function: functionFits adapted from https://www.statforbiology.com/nonlinearregression/usefulequations

Function description: fit raw agreement values to growth (for peakiness, agreement index, surprise index, detection accuracy and predictive value) or decay (for peak-to-peak distance) functions

Parameter y: array of dependent variable

Precondition y: numeric

Parameter x: array of independent variable

Precondition x: factor (note: this parameter is a factor because it was taken from data frame that’s fed into lmer)

Parameter sampSize: sequence of sample size

Precondition sampSize:

functionFit <- function(y, x, sampSize){
  x <- as.numeric(as.character(x))
  
  linear <- lm(y~x)
  exponential <- tryCatch(drm(y~x, fct = DRC.expoDecay()), error = function(e) drm(y~x, fct = EXD.2()))
  log <- drm(y ~ x, fct = DRC.logCurve())
  asympReg <- drm(y ~ x, fct = DRC.asymReg())
  powerCurve <- tryCatch(drm(y~x, fct = DRC.powerCurve()), error = function(e) drm(y~x, fct = L.4()))
  yl <- nls(y~NLS.YL(x, a, i)) #yielded erratic behavior for normativeHit
  
  models <- c('linear','exponential', 'log', 'asympReg', 'powerCurve', 'yl')
  model.prediction <- data.frame()
  for (model in models){
          assign("prediction", get(model))
          model.prediction <- rbind(model.prediction, data.frame(func = model, sampleSize = x, prediction = predict(prediction)))
  }
  
  bic <- BIC(linear,exponential,log,asympReg,powerCurve, yl)
  aic <- AIC(linear,exponential, log, asympReg, powerCurve, yl)
  
  functionFit.summary <- data.frame(AIC_df = aic$df,
                                    AIC = aic$AIC, 
                                    BIC_df = bic$df,
                                    BIC = bic$BIC)
  
  row.names(functionFit.summary) <- models
  
  mod.predicted <- predict(get(models[match(min(bic$BIC), bic$BIC)]), data.frame(sampSize))#, interval = 'confidence', level = .99
  elbow_pos <- elbow_finder(sampSize, mod.predicted)
  model <- models[match(min(bic$BIC), bic$BIC)]
  
  return(list(model.prediction, functionFit.summary,elbow_pos, model, coef(get(model))[3]))
}

Function: plotIndivAgreement

Function description: function to plot agreement values

Parameter df: data frame containing values for plotting

Precondition df: data frame

Parameter var: variable name for y values (i.e. which agreement measure)

Precondition var: string

Parameter ylab: y-axis label

Precondition ylab: string

Parameter ymin: minimum y axis value

Precondition ymin: numeric

Parameter ymax: maximum y axis value

Precondition ymax: numeric

Parameter title: plot title

Precondition title: string

plotIndivAgreement <- function(df, var, ylab, ymin, ymax, title){
  factors <- factor.order[(length(unique(df$testCond))%%3):3]
  factor.colors <- factor.color[(length(unique(df$testCond))%%3):3]
  factor.shapes <- factor.shape[(length(unique(df$testCond))%%3):3]
  factor.labels <- factor.label[(length(unique(df$testCond))%%3):3]
  transform(df, testCond=factor(testCond,levels=factors))
  elbow = mean(df$elbow)

    figure <- ggplot(df, aes(x = sampleSize, y = get(var, df), fill = testCond, color = testCond, shape = testCond))+
    geom_jitter(alpha = jitter.alpha, width = jitter.width)+
    stat_summary(geom = "errorbar", fun.data = "mean_cl_normal", width = errorbar.width, color = 'black', alpha = errorbar.alpha)+
    geom_point(stat = "summary", fun = "mean", size = point.size, alpha = point.alpha, color = 'black')+
    scale_fill_manual(values = factor.colors, labels = factor.labels)+
    scale_color_manual(values = factor.colors, labels = factor.labels)+
    scale_shape_manual(values = factor.shapes, labels = factor.labels)+
    #scale_x_continuous(breaks = seq(from = 2, to = 32, by = 2))+
    geom_text(color = "black", x = elbow/2, y =  ymin + (elbowPos.scale*(ymax - (ymin))), label = '^', size = elbow.size, show.legend = FALSE)+
    #facet_wrap(~grain, scales='fixed', labeller = labeller(grain = c("c" = "Coarse", "f" = "Fine")))+
    # geom_text(color = 'black', x = 6, y = -.2, label = "+", size = txt.size, show.legend = FALSE)+
    labs(y = ylab, x = "Sample Size", title = title)+
    ylim(ymin,ymax)+
    theme.esMethods
    #assign(paste('figure.', g, sep = ''), figure)`
    #ggsave(paste(paste('figure.', g, sep = ''), ".png", sep = ""))
    
    return(figure)
}

Function: plotFunctionFits

Function description: function to plot fitted growth/ decay functions

Parameter df: data frame containing values for plotting (i.e. function fit predictions)

Precondition df: data frame

Parameter original_df: data frame containing raw agremement values (to add as jitter)

Precondition df: data frame

Parameter var: original_df variable name for agreement values (i.e. which agreement measure)

Precondition var: string

Parameter ylab: y-axis label

Precondition ylab: string

Parameter ymin: minimum y axis value

Precondition ymin: numeric

Parameter ymax: maximum y axis value

plotFunctionFits <- function(df, original_df, var, ylab, ymin, ymax){
  figure <- ggplot(df,aes(x = sampleSize, y = prediction, color = func))+
  geom_jitter(data = original_df[original_df$testCond == "sample",], 
              aes(x = as.numeric(as.character(sampleSize)), y = get(var, original_df)), alpha = fit.jitter.alpha, width = fit.jitter.width, fill = dot.color, color = dot.color)+
  geom_line(size = fit.line.size, alpha = fit.line.alpha)+
  scale_color_manual(values = line.colors,
                     labels = c("Linear", "Exponential", "Logarithmic", "Asymptotic", "Power Curve", "Yield Loss"))+
  geom_point(data = original_df[original_df$testCond == "sample",], aes (x = as.numeric(as.character(sampleSize)), y = get(var, original_df)), stat = "summary", fun = "mean", shape = 21, size = 4, color = 'black', fill = dot.color)+
  scale_x_continuous(breaks = c(seq(from =2, to = 32, by = 2)))+
  ylim(ymin,ymax)+
  labs(y = ylab, x = "Sample size")+
  #facet_wrap(~grain, scales='fixed', labeller = labeller(grain = c("c" = "Coarse", "f" = "Fine")))+
  theme.esMethods
  
  return(figure)
}

Analyse data

Peakiness

Load dataset and remove outliers for evAct peakiness (value > 20*sd)

setwd('../data/bootstrapped')
peakiness.comm <- read.delim("esMethods_Peakiness_adjust_c0.1_f0.05_100Iterations_commercial.txt", head = TRUE)
peakiness.comm$sampleSize <- as.factor(peakiness.comm$sampleSize)
peakiness.comm$peakiness <- peakiness.comm$act.peakiness/peakiness.comm$min.peakiness
peakiness.comm$logPeakiness <- log(peakiness.comm$peakiness, base = 10)

peakiness.evAct <- read.delim("esMethods_Peakiness_adjust_c0.1_f0.05_100Iterations_evAct.txt", head = TRUE)
peakiness.evAct$peakiness <- peakiness.evAct$act.peakiness/peakiness.evAct$min.peakiness
peakiness.evAct <- peakiness.evAct[peakiness.evAct$peakiness < (mean(peakiness.evAct$peakiness) + (20*sd(peakiness.evAct$peakiness))),]
peakiness.evAct$sampleSize <- as.factor(peakiness.evAct$sampleSize)
peakiness.evAct$logActPeakiness <- log(peakiness.evAct$act.peakiness, base = 10)
peakiness.evAct$logPeakiness <- log(peakiness.evAct$peakiness, base = 10)

Fit peakiness values to growth and decay functions

##commercial 
peakiness.fit.coarse.comm <- functionFit(peakiness.comm$peakiness[peakiness.comm$grain == "c" & peakiness.comm$testCond == 'sample'], peakiness.comm$sampleSize[peakiness.comm$grain == "c" & peakiness.comm$testCond == 'sample'], sampleSize)
peakiness.fit.fine.comm <- functionFit(peakiness.comm$peakiness[peakiness.comm$grain == "f" & peakiness.comm$testCond == 'sample'], peakiness.comm$sampleSize[peakiness.comm$grain == "f" & peakiness.comm$testCond == 'sample'], sampleSize)

peakiness.functionFit <- peakiness.fit.coarse.comm[[1]] %>% mutate(grain = "c")
peakiness.functionFit <- rbind(peakiness.functionFit, peakiness.fit.fine.comm[[1]] %>% mutate(grain = "f"))
peakiness.functionFit$grain <- as.factor(peakiness.functionFit$grain)

#peakiness.functionFit.comm.figure <- plotFunctionFits(peakiness.functionFit, peakiness.comm[peakiness.comm$testCond == "sample",], "peakiness", "Peakiness (raw)")

#add elbow value to main dataframe 
for (i in 1:nrow(peakiness.comm)){
  if(peakiness.comm$grain[i] == 'c'){
    peakiness.comm$elbow[i] = 0 #coarse commercial peakiness did not stabilize
  } else {
    peakiness.comm$elbow[i] = peakiness.functionFit[[3]][1]
  }
}

##everyday Activities 
peakiness.fit.coarse.evAct <- functionFit(peakiness.evAct$peakiness[peakiness.evAct$grain == "c" & peakiness.evAct$testCond == "sample"], peakiness.evAct$sampleSize[peakiness.evAct$grain == 'c' & peakiness.evAct$testCond == "sample"], sampleSize)
peakiness.fit.fine.evAct <- functionFit(peakiness.evAct$peakiness[peakiness.evAct$grain == "f" & peakiness.evAct$testCond == "sample"], peakiness.evAct$sampleSize[peakiness.evAct$grain == 'f' & peakiness.evAct$testCond == "sample"], sampleSize)

peakiness.evAct.functionFit <- peakiness.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
peakiness.evAct.functionFit <- rbind(peakiness.evAct.functionFit, peakiness.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
peakiness.evAct.functionFit$grain <- as.factor(peakiness.evAct.functionFit$grain)

#peakiness.functionFit.evAct.figure <- plotFunctionFits(peakiness.evAct.functionFit, peakiness.evAct[peakiness.evAct$testCond == "sample",], "peakiness", "Peakiness (raw)") + ylim(0,50)

#add elbow value to main dataframe 
for (i in 1:nrow(peakiness.evAct)){
  if(peakiness.evAct$grain[i] == 'c'){
    peakiness.evAct$elbow[i] = peakiness.fit.coarse.evAct[[3]][1]
  } else {
    peakiness.evAct$elbow[i] = 0 #fine evAct peakiness did not stabilize
  }
}

Plot peakiness over sample size

peakiness.comm.coarse.plot <- plotIndivAgreement(peakiness.comm[peakiness.comm$grain == 'c',], "logPeakiness", "Peakiness (log)", 0, 2.5, "Coarse")
peakiness.comm.fine.plot <- plotIndivAgreement(peakiness.comm[peakiness.comm$grain == 'f',],"logPeakiness", "Peakiness (log)", 0, 2.5, "Fine")

peakiness.evAct.coarse.plot <- plotIndivAgreement(peakiness.evAct[peakiness.evAct$grain == 'c', ], "logPeakiness", "Peakiness (log)", 0, 2.5, "")
peakiness.evAct.fine.plot <- plotIndivAgreement(peakiness.evAct[peakiness.evAct$grain == 'f', ],  "logPeakiness", "Peakiness (log)", 0, 2.5, "")

peakiness.plot.combined <- ggarrange(peakiness.comm.coarse.plot, peakiness.comm.fine.plot,
                                      peakiness.evAct.coarse.plot, peakiness.evAct.fine.plot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")
peakiness.plot.combined

#ggsave("esMethod_peakiness.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")
#ggsave("esMethod_peakiness.eps", device = cairo_ps, width = 25, height = 19, units = "cm", bg = "transparent")

Plot peakiness function fits

peakiness.comm.coarse.funcPlot <- plotFunctionFits(peakiness.functionFit[peakiness.functionFit$grain == 'c',], peakiness.comm[peakiness.comm$grain == 'c' & peakiness.comm$testCond == 'sample',], "peakiness", "Peakiness", 2, 24)
peakiness.comm.fine.funcPlot <- plotFunctionFits(peakiness.functionFit[peakiness.functionFit$grain == 'f',], peakiness.comm[peakiness.comm$grain == 'f' & peakiness.comm$testCond == 'sample',], "peakiness", "Peakiness", 2, 24)

peakiness.evAct.coarse.funcPlot <- plotFunctionFits(peakiness.evAct.functionFit[peakiness.evAct.functionFit$grain == 'c',], peakiness.evAct[peakiness.evAct$grain == 'c' & peakiness.evAct$testCond == 'sample',], "peakiness", "Peakiness", 2, 50)
peakiness.evAct.fine.funcPlot <- plotFunctionFits(peakiness.evAct.functionFit[peakiness.evAct.functionFit$grain == 'f',], peakiness.evAct[peakiness.evAct$grain == 'f' & peakiness.evAct$testCond == 'sample',], "peakiness", "Peakiness", 2, 50)

peakiness.funcFit.plot.combined <- ggarrange(peakiness.comm.coarse.funcPlot, peakiness.comm.fine.funcPlot,
                                      peakiness.evAct.coarse.funcPlot, peakiness.evAct.fine.funcPlot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

peakiness.funcFit.plot.combined

#ggsave("esMethod_peakiness_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")

Peakiness linear models and analysis

peakiness.comm.lmer <- lmer(peakiness~sampleSize*grain*testCond + (1|movieName), data = peakiness.comm)
#emm_options(pbkrtest.limit = 3000)
peakiness.comm.samp_vs_rand <- summary(emmeans(peakiness.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1))), infer = T)
peakiness.comm.pairwise <- summary(emmeans(peakiness.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

peakiness.comm.poly <- summary(contrast(emmeans(peakiness.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast

peakiness.evAct.lmer <- lmer(peakiness~sampleSize*grain*testCond + (1|movieName), data = peakiness.evAct)
#emm_options(lmerTest.limit = 25600)
peakiness.evAct.samp_vs_rand <- summary(emmeans(peakiness.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1))), infer = T)# overall random vs sample 
peakiness.evAct.pairwise <- summary(emmeans(peakiness.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#pairwise comparison 
peakiness.evAct.poly <- summary(contrast(emmeans(peakiness.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast

Peakiness Sample - random contrast commercial-lab

peakiness.comm.samp_vs_rand$contrasts
## grain = c:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       4.587464 0.03646846 Inf  4.515987  4.658941 125.793  <.0001
## 
## grain = f:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       5.859044 0.03646846 Inf  5.787567  5.930521 160.661  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peakiness sample-random contrast everyday-online

peakiness.evAct.samp_vs_rand$contrasts
## grain = c:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       15.657130 0.2507069 Inf 15.165753 16.148506  62.452  <.0001
## 
## grain = f:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1        6.619238 0.2506576 Inf  6.127958  7.110518  26.407  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peakiness smallest pairwise effect commercial-lab

peakiness.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 2 x 10
## # Groups:   grain [2]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 6          c        1.34  0.146   Inf    1.05       1.63 
## 2 sample … 4          f        0.347 0.146   Inf    0.0614     0.633
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Peakiness smallest pairwise effect everyday-online

peakiness.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 2 x 10
## # Groups:   grain [2]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c        10.7   1.00   Inf      8.75     12.7 
## 2 sample … 2          f         3.68  1.00   Inf      1.72      5.65
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Peakiness polynomial fit commercial-lab

peakiness.comm.poly[peakiness.comm.poly$testCond == "sample",]
##     contrast grain testCond    estimate         SE  df  asymp.LCL
## 7     linear     c   sample  309.806809   3.803925 Inf  302.35125
## 8  quadratic     c   sample   33.310526   7.795726 Inf   18.03118
## 9      cubic     c   sample    4.709911 103.547837 Inf -198.24012
## 10    linear     f   sample  336.377790   3.803925 Inf  328.92223
## 11 quadratic     f   sample  -76.594747   7.795726 Inf  -91.87409
## 12     cubic     f   sample -302.832594 103.547837 Inf -505.78262
##    asymp.UCL     z.ratio      p.value
## 7  317.26236 81.44398802 0.000000e+00
## 8   48.58987  4.27292157 1.929283e-05
## 9  207.65994  0.04548536 9.637204e-01
## 10 343.83335 88.42913694 0.000000e+00
## 11 -61.31540 -9.82522292 8.767983e-23
## 12 -99.88256 -2.92456707 3.449357e-03

Peakiness polynomial fit commercial-lab

peakiness.evAct.poly[peakiness.evAct.poly$testCond == "sample",]
##     contrast grain testCond   estimate        SE  df  asymp.LCL  asymp.UCL
## 7     linear     c   sample  197.32892  26.16873 Inf  146.03915  248.61868
## 8  quadratic     c   sample -146.86548  53.63168 Inf -251.98164  -41.74933
## 9      cubic     c   sample 1707.65861 712.29000 Inf  311.59586 3103.72135
## 10    linear     f   sample  -80.43273  26.14541 Inf -131.67678  -29.18868
## 11 quadratic     f   sample   32.04620  53.58214 Inf  -72.97286  137.06526
## 12     cubic     f   sample  752.69029 711.71232 Inf -642.24022 2147.62081
##       z.ratio      p.value
## 7   7.5406381 4.676775e-14
## 8  -2.7384094 6.173718e-03
## 9   2.3974204 1.651097e-02
## 10 -3.0763620 2.095433e-03
## 11  0.5980762 5.497891e-01
## 12  1.0575766 2.902485e-01

Peak-to-peak distance

Load dataset

setwd('../data/bootstrapped')
peaktopeak.comm <- read.delim("esMethods_PeakPeakDistance_c0.1_f0.05_100Iterations_commercial.txt", head = TRUE)
peaktopeak.comm$sampleSize <- as.factor(peaktopeak.comm$sampleSize)
peaktopeak.comm$logPeaktopeak <- log(peaktopeak.comm$nearestDistance, base = 10)

peaktopeak.evAct <- read.delim("esMethods_PeakPeakDistance_c0.1_f0.05_100Iterations_evAct.txt", head = TRUE)
peaktopeak.evAct$sampleSize <- as.factor(peaktopeak.evAct$sampleSize)
peaktopeak.evAct$logPeaktopeak <- log(peaktopeak.evAct$nearestDistance, base = 10)

Fit peak-to-peak values to growth and decay functions

##commercial 
peaktopeak.fit.coarse.comm <- functionFit(peaktopeak.comm$nearestDistance[peaktopeak.comm$grain == "c" & peaktopeak.comm$testCond == 'sample'], peaktopeak.comm$sampleSize[peaktopeak.comm$grain == "c" & peaktopeak.comm$testCond == 'sample'], sampleSize)
peaktopeak.fit.fine.comm <- functionFit(peaktopeak.comm$nearestDistance[peaktopeak.comm$grain == "f" & peaktopeak.comm$testCond == 'sample'], peaktopeak.comm$sampleSize[peaktopeak.comm$grain == "f" & peaktopeak.comm$testCond == 'sample'], sampleSize)

peaktopeak.functionFit <- peaktopeak.fit.coarse.comm[[1]] %>% mutate(grain = "c")
peaktopeak.functionFit <- rbind(peaktopeak.functionFit, peaktopeak.fit.fine.comm[[1]] %>% mutate(grain = "f"))
peaktopeak.functionFit$grain <- as.factor(peaktopeak.functionFit$grain)

#peaktopeak.functionFit.comm.figure <- plotFunctionFits(peaktopeak.functionFit, peaktopeak.comm[peaktopeak.comm$testCond == "sample",], "nearestDistance", "Peak-to-peak distance (s)")

for(i in 1:nrow(peaktopeak.comm)){
  if(peaktopeak.comm$grain[i] == 'c'){
    peaktopeak.comm$elbow[i] = peaktopeak.fit.coarse.comm[[3]][1]
  }else{
    peaktopeak.comm$elbow[i] = peaktopeak.fit.fine.comm[[3]][1]
  }
}

##everyday Activities 
peaktopeak.fit.coarse.evAct <- functionFit(peaktopeak.evAct$nearestDistance[peaktopeak.evAct$grain == "c" & peaktopeak.evAct$testCond == "sample"], peaktopeak.evAct$sampleSize[peaktopeak.evAct$grain == 'c' & peaktopeak.evAct$testCond == "sample"], sampleSize)
peaktopeak.fit.fine.evAct <- functionFit(peaktopeak.evAct$nearestDistance[peaktopeak.evAct$grain == "f" & peaktopeak.evAct$testCond == "sample"], peaktopeak.evAct$sampleSize[peaktopeak.evAct$grain == 'f' & peaktopeak.evAct$testCond == "sample"], sampleSize)

peaktopeak.evAct.functionFit <- peaktopeak.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
peaktopeak.evAct.functionFit <- rbind(peaktopeak.evAct.functionFit, peaktopeak.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
peaktopeak.evAct.functionFit$grain <- as.factor(peaktopeak.evAct.functionFit$grain)

#peaktopeak.functionFit.evAct.figure <- plotFunctionFits(peaktopeak.evAct.functionFit, peaktopeak.evAct[peaktopeak.evAct$testCond == "sample",], "nearestDistance", "Peak-to-peak distance (s)") + ylim(0,50)

for(i in 1:nrow(peaktopeak.evAct)){
  if(peaktopeak.evAct$grain[i] == 'c'){
    peaktopeak.evAct$elbow[i] = peaktopeak.fit.coarse.evAct[[3]][1]
  }else{
    peaktopeak.evAct$elbow[i] = peaktopeak.fit.fine.evAct[[3]][1]
  }
}

Plot peak-to-peak over sample size

Log10 transformed peak-to-peak distance over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

peaktopeak.comm.coarse.plot <- plotIndivAgreement(peaktopeak.comm[peaktopeak.comm$grain == 'c',], "logPeaktopeak", "Peak to peak distance (log)", 0, 2, "Coarse")
peaktopeak.comm.fine.plot <- plotIndivAgreement(peaktopeak.comm[peaktopeak.comm$grain == 'f',], "logPeaktopeak", "Peak to peak distance (log)", -0.5, 1.5, "Fine")

peaktopeak.evAct.coarse.plot <- plotIndivAgreement(peaktopeak.evAct[peaktopeak.evAct$grain == 'c', ], "logPeaktopeak", "Peak to peak distance (log)", 0, 2, "")
peaktopeak.evAct.fine.plot <- plotIndivAgreement(peaktopeak.evAct[peaktopeak.evAct$grain == 'f', ], "logPeaktopeak", "Peak to peak distance (log)", -0.5, 1.5, "")

peaktopeak.plot.combined <- ggarrange(peaktopeak.comm.coarse.plot, peaktopeak.comm.fine.plot,
                                      peaktopeak.evAct.coarse.plot, peaktopeak.evAct.fine.plot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

peaktopeak.plot.combined

setwd('../plots')
ggsave("esMethod_peaktopeak.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent", dpi = 300)
# ggsave(plot = peaktopeak.plot.combined, file = "esMethod_peaktopeak.png", 
#        type = "cairo-png",  bg = "transparent",
#        width = 25, height = 19, units = "cm", dpi = 300)

Plot peak-to-peak function fits

Peak-to-peak values fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

#df, original_df, var, ylab, ymin, ymax, asymptote
peaktopeak.comm.coarse.funcPlot <- plotFunctionFits(peaktopeak.functionFit[peaktopeak.functionFit$grain == 'c',], peaktopeak.comm[peaktopeak.comm$grain == 'c' & peaktopeak.comm$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 15)
peaktopeak.comm.fine.funcPlot <- plotFunctionFits(peaktopeak.functionFit[peaktopeak.functionFit$grain == 'f',], peaktopeak.comm[peaktopeak.comm$grain == 'f' & peaktopeak.comm$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 10)

peaktopeak.evAct.coarse.funcPlot <- plotFunctionFits(peaktopeak.evAct.functionFit[peaktopeak.evAct.functionFit$grain == 'c',], peaktopeak.evAct[peaktopeak.evAct$grain == 'c' & peaktopeak.evAct$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 25)
peaktopeak.evAct.fine.funcPlot <- plotFunctionFits(peaktopeak.evAct.functionFit[peaktopeak.evAct.functionFit$grain == 'f',], peaktopeak.evAct[peaktopeak.evAct$grain == 'f' & peaktopeak.evAct$testCond == 'sample',], "nearestDistance", "Peak-to-peak distance (s)", 0, 5)

peaktopeak.funcFit.plot.combined <- ggarrange(peaktopeak.comm.coarse.funcPlot, peaktopeak.comm.fine.funcPlot,
                                      peaktopeak.evAct.coarse.funcPlot, peaktopeak.evAct.fine.funcPlot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

peaktopeak.funcFit.plot.combined

setwd('../plots')
ggsave("esMethod_peaktopeak_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")

Peak-to-peak linear models and analysis

#commercial-lab 
peaktopeak.comm.lmer <- lmer(nearestDistance~sampleSize*grain*testCond + (1|movieName), data = peaktopeak.comm)
#emm_options(pbkrtest.limit = 19200)
peaktopeak.comm.poly <- summary(contrast(emmeans(peaktopeak.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
peaktopeak.comm.rand_vs_samp <- summary(emmeans(peaktopeak.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, 1, -1))), infer = T)# overall random - sample 
peaktopeak.comm.crossMov_vs_samp <- summary(emmeans(peaktopeak.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, 0, -1))), infer = T)# overall crossMov - sample 
peaktopeak.comm.crossMov_vs_rand <- summary(emmeans(peaktopeak.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random 

peaktopeak.comm.pairwise <- summary(emmeans(peaktopeak.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "pairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

#everyday-online
peaktopeak.evAct.lmer <- lmer(nearestDistance~sampleSize*grain*testCond + (1|movieName), data = peaktopeak.evAct)
#emm_options(pbkrtest.limit = 3000)
peaktopeak.evAct.poly <- summary(contrast(emmeans(peaktopeak.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
peaktopeak.evAct.rand_vs_samp <- summary(emmeans(peaktopeak.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, 1, -1))), infer = T)# overall random - sample 
peaktopeak.evAct.crossMov_vs_samp <- summary(emmeans(peaktopeak.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, 0, -1))), infer = T)# overall crossMov - sample 
peaktopeak.evAct.crossMov_vs_rand <- summary(emmeans(peaktopeak.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random 

peaktopeak.evAct.pairwise <- summary(emmeans(peaktopeak.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "pairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

Peak-to-peak random - sample contrast commercial-lab

peaktopeak.comm.rand_vs_samp$contrasts
## grain = c:
##  contrast  estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       3.0491790 0.03624087 Inf 2.9781482 3.1202098  84.136  <.0001
## 
## grain = f:
##  contrast  estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.7547491 0.03624087 Inf 0.6837183 0.8257799  20.826  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peak-to-peak random - sample contrast everyday-online

peaktopeak.evAct.rand_vs_samp$contrasts
## grain = c:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       2.6898060 0.1011654 Inf 2.4915256 2.8880865  26.588  <.0001
## 
## grain = f:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.3761267 0.1011654 Inf 0.1778462 0.5744071   3.718  0.0002
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peak-to-peak cross movie - sample contrast commercial-lab

peaktopeak.comm.crossMov_vs_samp$contrasts
## grain = c:
##  contrast  estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       3.3981889 0.03624087 Inf  3.327158 3.4692197  93.767  <.0001
## 
## grain = f:
##  contrast  estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.7305018 0.03624087 Inf  0.659471 0.8015326  20.157  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peak-to-peak cross movie - sample contrast everyday-online

peaktopeak.evAct.crossMov_vs_samp$contrasts
## grain = c:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       5.3137203 0.1011654 Inf 5.1154399 5.5120008  52.525  <.0001
## 
## grain = f:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.7301477 0.1011654 Inf 0.5318673 0.9284282   7.217  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peak-to-peak cross movie - random contrast commercial-lab

peaktopeak.comm.crossMov_vs_rand$contrasts
## grain = c:
##  contrast    estimate         SE  df   asymp.LCL  asymp.UCL z.ratio
##  c1        0.34900988 0.03624087 Inf  0.27797908 0.42004068   9.630
##  p.value
##   <.0001
## 
## grain = f:
##  contrast    estimate         SE  df   asymp.LCL  asymp.UCL z.ratio
##  c1       -0.02424736 0.03624087 Inf -0.09527815 0.04678344  -0.669
##  p.value
##   0.5035
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peak-to-peak cross movie - random contrast everyday-online

peaktopeak.evAct.crossMov_vs_rand$contrasts
## grain = c:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       2.6239143 0.1011654 Inf 2.4256338 2.8221948  25.937  <.0001
## 
## grain = f:
##  contrast  estimate        SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.3540211 0.1011654 Inf 0.1557406 0.5523015   3.499  0.0005
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Peak-to-peak smallest sample size for significant pairwise effect commercial-lab

peaktopeak.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 5 x 10
## # Groups:   grain, contrast [5]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 crossMo… 2          c        0.295 0.145   Inf    0.0111     0.579
## 2 crossMo… 2          c        2.70  0.145   Inf    2.42       2.98 
## 3 random … 2          c        2.40  0.145   Inf    2.12       2.69 
## 4 crossMo… 2          f        0.557 0.145   Inf    0.273      0.841
## 5 random … 2          f        0.726 0.145   Inf    0.441      1.01 
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Peak-to-peak smallest sample size for significant pairwise effect everyday-online

peaktopeak.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 3 x 10
## # Groups:   grain, contrast [3]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 crossMo… 2          c         1.25 0.405   Inf     0.457      2.04
## 2 crossMo… 2          c         3.79 0.405   Inf     3.00       4.58
## 3 random … 2          c         2.54 0.405   Inf     1.75       3.33
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Peak-to-peak smallest pairwise effect commercial-lab

peaktopeak.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 5 x 10
## # Groups:   grain, contrast [5]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 crossMo… 2          c        0.295 0.145   Inf    0.0111     0.579
## 2 crossMo… 2          c        2.70  0.145   Inf    2.42       2.98 
## 3 random … 4          c        2.10  0.145   Inf    1.82       2.38 
## 4 crossMo… 2          f        0.557 0.145   Inf    0.273      0.841
## 5 random … 12         f        0.642 0.145   Inf    0.358      0.927
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Peak-to-peak smallest pairwise effect everyday-online coarse

peaktopeak.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 3 x 10
## # Groups:   grain, contrast [3]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 crossMo… 2          c         1.25 0.405   Inf     0.457      2.04
## 2 crossMo… 2          c         3.79 0.405   Inf     3.00       4.58
## 3 random … 4          c         2.09 0.405   Inf     1.30       2.88
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Peak-to-peak largest pairwise effect everyday-online fine (note: non significant effect)

peaktopeak.evAct.pairwise %>% dplyr::filter(p.value > .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == max(z.ratio))
## # A tibble: 3 x 10
## # Groups:   grain, contrast [3]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 random … 22         f        0.433 0.405   Inf  -3.60e-1      1.23
## 2 crossMo… 32         f        0.408 0.405   Inf  -3.85e-1      1.20
## 3 crossMo… 32         f        0.793 0.405   Inf  -5.30e-5      1.59
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Peak-to-peak polynomial fit commercial-lab

peaktopeak.comm.poly[peaktopeak.comm.poly$testCond == "sample",]
##     contrast grain testCond    estimate         SE  df   asymp.LCL
## 13    linear     c   sample  -215.37582   3.780185 Inf  -222.78484
## 14 quadratic     c   sample   295.65601   7.747074 Inf   280.47203
## 15     cubic     c   sample -2278.86403 102.901608 Inf -2480.54747
## 16    linear     f   sample   -66.42806   3.780185 Inf   -73.83709
## 17 quadratic     f   sample    89.66520   7.747074 Inf    74.48122
## 18     cubic     f   sample  -711.27868 102.901608 Inf  -912.96212
##      asymp.UCL    z.ratio       p.value
## 13  -207.96679 -56.974941  0.000000e+00
## 14   310.84000  38.163572  0.000000e+00
## 15 -2077.18058 -22.146049 1.138796e-108
## 16   -59.01904 -17.572702  3.987230e-69
## 17   104.84919  11.574073  5.577059e-31
## 18  -509.59523  -6.912221  4.771230e-12

Peak-to-peak polynomial fit everyday-online

peaktopeak.evAct.poly[peaktopeak.evAct.poly$testCond == "sample",]
##     contrast grain testCond    estimate        SE  df   asymp.LCL
## 13    linear     c   sample  -144.34337  10.55228 Inf  -165.02546
## 14 quadratic     c   sample   200.56105  21.62574 Inf   158.17538
## 15     cubic     c   sample -2036.14458 287.24691 Inf -2599.13818
## 16    linear     f   sample   -32.64406  10.55228 Inf   -53.32615
## 17 quadratic     f   sample    53.53322  21.62574 Inf    11.14756
## 18     cubic     f   sample  -451.80064 287.24691 Inf -1014.79424
##      asymp.UCL    z.ratio      p.value
## 13  -123.66128 -13.678881 1.357655e-42
## 14   242.94671   9.274183 1.789848e-20
## 15 -1473.15098  -7.088482 1.355909e-12
## 16   -11.96197  -3.093556 1.977736e-03
## 17    95.91889   2.475440 1.330719e-02
## 18   111.19296  -1.572865 1.157501e-01

Agreement Index

Load dataset

setwd('../data/bootstrapped')
agreement.comm <- read.delim("esMethods_AgreementIndex_100Iterations_commercial.txt", head = TRUE)
agreement.comm$sampleSize <- as.factor(agreement.comm$sampleSize)

agreement.evAct <- read.delim("esMethods_AgreementIndex_100Iterations_evAct.txt", head = TRUE)
agreement.evAct$sampleSize <- as.factor(agreement.evAct$sampleSize)

Fit Agreement Index values to growth and decay functions

##commercial-lab 
agreement.fit.coarse.comm <- functionFit(agreement.comm$agreementIndex[agreement.comm$grain == "c" & agreement.comm$testCond == 'sample'], agreement.comm$sampleSize[agreement.comm$grain == "c" & agreement.comm$testCond == 'sample'], sampleSize)
agreement.fit.fine.comm <- functionFit(agreement.comm$agreementIndex[agreement.comm$grain == "f" & agreement.comm$testCond == 'sample'], agreement.comm$sampleSize[agreement.comm$grain == "f" & agreement.comm$testCond == 'sample'], sampleSize)

agreement.functionFit <- agreement.fit.coarse.comm[[1]] %>% mutate(grain = "c")
agreement.functionFit <- rbind(agreement.functionFit, agreement.fit.fine.comm[[1]] %>% mutate(grain = "f"))
agreement.functionFit$grain <- as.factor(agreement.functionFit$grain)

for(i in 1:nrow(agreement.comm)){
  if(agreement.comm$grain[i] == 'c'){
    agreement.comm$elbow[i] = agreement.fit.coarse.comm[[3]][1]
  } else {
    agreement.comm$elbow[i] = agreement.fit.fine.comm[[3]][1]
  }
}

##everyday Activities-lab 
agreement.fit.coarse.evAct <- functionFit(agreement.evAct$agreementIndex[agreement.evAct$grain == "c" & agreement.evAct$testCond == "sample"], agreement.evAct$sampleSize[agreement.evAct$grain == 'c' & agreement.evAct$testCond == "sample"], sampleSize)
agreement.fit.fine.evAct <- functionFit(agreement.evAct$agreementIndex[agreement.evAct$grain == "f" & agreement.evAct$testCond == "sample"], agreement.evAct$sampleSize[agreement.evAct$grain == 'f' & agreement.evAct$testCond == "sample"], sampleSize)

agreement.evAct.functionFit <- agreement.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
agreement.evAct.functionFit <- rbind(agreement.evAct.functionFit, agreement.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
agreement.evAct.functionFit$grain <- as.factor(agreement.evAct.functionFit$grain)

for(i in 1:nrow(agreement.evAct)){
  if(agreement.evAct$grain[i] == 'c'){
    agreement.evAct$elbow[i] = agreement.fit.coarse.evAct[[3]][1]
  } else {
    agreement.evAct$elbow[i] = agreement.fit.fine.evAct[[3]][1]
  }
}

Plot Agreement Index over sample size

Agreement Index over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

agreement.comm.coarse.plot <- plotIndivAgreement(agreement.comm[agreement.comm$grain == 'c',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1, "Coarse")
agreement.comm.fine.plot <- plotIndivAgreement(agreement.comm[agreement.comm$grain == 'f',],"agreementIndex", "Agreement Index", -0.0000000000000001, 1, "Fine")

agreement.evAct.coarse.plot <- plotIndivAgreement(agreement.evAct[agreement.evAct$grain == 'c', ], "agreementIndex", "Agreement Index", -0.0000000000000001, 1, "")
agreement.evAct.fine.plot <- plotIndivAgreement(agreement.evAct[agreement.evAct$grain == 'f', ],  "agreementIndex", "Agreement Index", -0.0000000000000001, 1, "")

agreement.plot.combined <- ggarrange(agreement.comm.coarse.plot, agreement.comm.fine.plot,
                                      agreement.evAct.coarse.plot, agreement.evAct.fine.plot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

agreement.plot.combined

setwd('../plots')
ggsave("esMethod_agreement.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")
# ggsave(plot = agreement.plot.combined, file = "esMethod_agreement.png", 
#        type = "cairo-png",  bg = "transparent",
#        width = 25, height = 19, units = "cm", dpi = 300)

Plot Agreement Index function fits

Agreement Index fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

agreement.comm.coarse.funcPlot <- plotFunctionFits(agreement.functionFit[agreement.functionFit$grain == 'c',], agreement.comm[agreement.comm$grain == 'c' & agreement.comm$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)
agreement.comm.fine.funcPlot <- plotFunctionFits(agreement.functionFit[agreement.functionFit$grain == 'f',], agreement.comm[agreement.comm$grain == 'f' & agreement.comm$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)

agreement.evAct.coarse.funcPlot <- plotFunctionFits(agreement.evAct.functionFit[agreement.evAct.functionFit$grain == 'c',], agreement.evAct[agreement.evAct$grain == 'c' & agreement.evAct$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)
agreement.evAct.fine.funcPlot <- plotFunctionFits(agreement.evAct.functionFit[agreement.evAct.functionFit$grain == 'f',], agreement.evAct[agreement.evAct$grain == 'f' & agreement.evAct$testCond == 'sample',], "agreementIndex", "Agreement Index", -0.0000000000000001, 1)

agreement.funcFit.plot.combined <- ggarrange(agreement.comm.coarse.funcPlot, agreement.comm.fine.funcPlot,
                                      agreement.evAct.coarse.funcPlot, agreement.evAct.fine.funcPlot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

agreement.funcFit.plot.combined

setwd('../plots')
ggsave("esMethod_agreement_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")

Agreement Index linear models and analysis

#Commercial 
agreement.comm.lmer <- lmer(agreementIndex~sampleSize*grain*testCond + (1|movieName), data = agreement.comm)
#emm_options(pbkrtest.limit = 3000)
agreement.comm.poly <- summary(contrast(emmeans(agreement.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
agreement.comm.rand_vs_samp <- summary(emmeans(agreement.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample 
agreement.comm.crossMov_vs_samp <- summary(emmeans(agreement.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample 
agreement.comm.crossMov_vs_rand <- summary(emmeans(agreement.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1, 0))), infer = T)# overall crossMov - random 

agreement.comm.pairwise <- summary(emmeans(agreement.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

agreement.evAct.lmer <- lmer(agreementIndex~sampleSize*grain*testCond + (1|movieName), data = agreement.evAct)
#emm_options(lmerTest.limit = 38688)
agreement.evAct.poly <- summary(contrast(emmeans(agreement.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
agreement.evAct.rand_vs_samp <- summary(emmeans(agreement.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample 
agreement.evAct.crossMov_vs_samp <- summary(emmeans(agreement.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample 
agreement.evAct.crossMov_vs_rand <- summary(emmeans(agreement.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1, 0))), infer = T)# overall crossMov - random 

agreement.evAct.pairwise <- summary(emmeans(agreement.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

Agreement Index random - sample contrast commercial-lab

agreement.comm.rand_vs_samp$contrasts
## grain = c:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.1577242 0.001170359 Inf 0.1554303 0.1600180 134.766  <.0001
## 
## grain = f:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.1205933 0.001170359 Inf 0.1182995 0.1228872 103.040  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Agreement Index random - sample contrast everyday-online

agreement.evAct.rand_vs_samp$contrasts
## grain = c:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.1904387 0.001051004 Inf 0.1883788 0.1924986 181.197  <.0001
## 
## grain = f:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.1612430 0.001051004 Inf 0.1591831 0.1633029 153.418  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Agreement Index cross movie - sample contrast commercial-lab

agreement.comm.crossMov_vs_samp$contrasts
## grain = c:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.2291428 0.001170359 Inf 0.2268489 0.2314366 195.788  <.0001
## 
## grain = f:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.2037672 0.001170359 Inf 0.2014733 0.2060611 174.107  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Agreement Index cross movie - sample contrast everyday-online

agreement.evAct.crossMov_vs_samp$contrasts
## grain = c:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.2232330 0.001051004 Inf 0.2211731 0.2252929  212.40  <.0001
## 
## grain = f:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.1864058 0.001051004 Inf 0.1843458 0.1884657  177.36  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Agreement Index cross movie - random contrast commercial-lab

agreement.comm.crossMov_vs_rand$contrasts
## grain = c:
##  contrast   estimate          SE  df  asymp.LCL  asymp.UCL z.ratio p.value
##  c1       0.07141858 0.001170359 Inf 0.06912472 0.07371244  61.023  <.0001
## 
## grain = f:
##  contrast   estimate          SE  df  asymp.LCL  asymp.UCL z.ratio p.value
##  c1       0.08317385 0.001170359 Inf 0.08087999 0.08546771  71.067  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Agreement Index cross movie - random contrast everyday-online

agreement.evAct.crossMov_vs_rand$contrasts
## grain = c:
##  contrast   estimate          SE  df  asymp.LCL  asymp.UCL z.ratio p.value
##  c1       0.03279429 0.001051004 Inf 0.03073436 0.03485422  31.203  <.0001
## 
## grain = f:
##  contrast   estimate          SE  df  asymp.LCL  asymp.UCL z.ratio p.value
##  c1       0.02516278 0.001051004 Inf 0.02310285 0.02722271  23.942  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Agreement Index smallest sample size for significant pairwise effect commercial-lab

agreement.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 6 x 10
## # Groups:   grain, contrast [6]
##   contrast sampleSize grain estimate      SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>   <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.0566 0.00468   Inf   0.0474     0.0658
## 2 sample … 2          c       0.0898 0.00468   Inf   0.0806     0.0989
## 3 random … 6          c       0.0129 0.00468   Inf   0.00370    0.0220
## 4 sample … 2          f       0.0400 0.00468   Inf   0.0308     0.0491
## 5 sample … 2          f       0.101  0.00468   Inf   0.0918     0.110 
## 6 random … 4          f       0.0117 0.00468   Inf   0.00248    0.0208
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Agreement Index smallest sample size for significant pairwise effect everyday-online

agreement.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 6 x 10
## # Groups:   grain, contrast [6]
##   contrast sampleSize grain estimate      SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>   <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.0469 0.00420   Inf   0.0386     0.0551
## 2 sample … 2          c       0.0803 0.00420   Inf   0.0721     0.0886
## 3 random … 12         c       0.0169 0.00420   Inf   0.00862    0.0251
## 4 sample … 2          f       0.0659 0.00420   Inf   0.0576     0.0741
## 5 sample … 2          f       0.134  0.00420   Inf   0.126      0.142 
## 6 random … 6          f       0.0109 0.00420   Inf   0.00262    0.0191
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Agreement Index smallest pairwise effect commercial-lab

agreement.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 6 x 10
## # Groups:   grain, contrast [6]
##   contrast sampleSize grain estimate      SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>   <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.0566 0.00468   Inf   0.0474     0.0658
## 2 sample … 2          c       0.0898 0.00468   Inf   0.0806     0.0989
## 3 random … 6          c       0.0129 0.00468   Inf   0.00370    0.0220
## 4 sample … 2          f       0.0400 0.00468   Inf   0.0308     0.0491
## 5 sample … 2          f       0.101  0.00468   Inf   0.0918     0.110 
## 6 random … 4          f       0.0117 0.00468   Inf   0.00248    0.0208
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Agreement Index smallest pairwise effect everyday-online

agreement.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 6 x 10
## # Groups:   grain, contrast [6]
##   contrast sampleSize grain estimate      SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>   <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.0469 0.00420   Inf   0.0386     0.0551
## 2 sample … 2          c       0.0803 0.00420   Inf   0.0721     0.0886
## 3 random … 12         c       0.0169 0.00420   Inf   0.00862    0.0251
## 4 sample … 2          f       0.0659 0.00420   Inf   0.0576     0.0741
## 5 sample … 2          f       0.134  0.00420   Inf   0.126      0.142 
## 6 random … 6          f       0.0109 0.00420   Inf   0.00262    0.0191
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Agreement Index polynomial fit commercial-lab

agreement.comm.poly[agreement.comm.poly$testCond == "sample",]
##     contrast grain testCond  estimate        SE  df asymp.LCL  asymp.UCL
## 13    linear     c   sample  10.18298 0.1220769 Inf   9.94371  10.422243
## 14 quadratic     c   sample -10.15537 0.2501832 Inf -10.64572  -9.665016
## 15     cubic     c   sample  71.02177 3.3230932 Inf  64.50863  77.534917
## 16    linear     f   sample  10.38112 0.1220769 Inf  10.14185  10.620384
## 17 quadratic     f   sample -11.83477 0.2501832 Inf -12.32512 -11.344421
## 18     cubic     f   sample  90.35719 3.3230932 Inf  83.84405  96.870336
##      z.ratio       p.value
## 13  83.41446  0.000000e+00
## 14 -40.59173  0.000000e+00
## 15  21.37219 2.424961e-101
## 16  85.03754  0.000000e+00
## 17 -47.30443  0.000000e+00
## 18  27.19069 8.369360e-163

Agreement Index polynomial fit everyday-online

agreement.evAct.poly[agreement.evAct.poly$testCond == "sample",]
##     contrast grain testCond   estimate        SE  df  asymp.LCL  asymp.UCL
## 13    linear     c   sample  11.456484 0.1096273 Inf  11.241618  11.671350
## 14 quadratic     c   sample -10.784017 0.2246692 Inf -11.224360 -10.343673
## 15     cubic     c   sample  70.339801 2.9842000 Inf  64.490877  76.188726
## 16    linear     f   sample   8.101342 0.1096273 Inf   7.886476   8.316208
## 17 quadratic     f   sample  -9.983328 0.2246692 Inf -10.423671  -9.542984
## 18     cubic     f   sample  78.066083 2.9842000 Inf  72.217158  83.915007
##      z.ratio       p.value
## 13 104.50391  0.000000e+00
## 14 -47.99954  0.000000e+00
## 15  23.57074 7.694444e-123
## 16  73.89893  0.000000e+00
## 17 -44.43568  0.000000e+00
## 18  26.15980 7.623629e-151

Surprise Index

Load dataset

setwd('../data/bootstrapped')

surprise.comm <- read.delim("esMethods_SurpriseIndex_c0.1_f0.05_100Iterations_commercial.txt", head = TRUE)
surprise.comm$sampleSize <- as.factor(surprise.comm$sampleSize)
surprise.comm$logSurprise <- log(surprise.comm$surpriseIndex, base = 10)

surprise.evAct <- read.delim("esMethods_SurpriseIndex_c0.1_f0.05_100Iterations_evAct.txt", head = TRUE)
surprise.evAct$sampleSize <- as.factor(surprise.evAct$sampleSize)
surprise.evAct$logSurprise <- log(surprise.evAct$surpriseIndex, base = 10)

Fit Surprise Index values to growth and decay functions

##commercial 
surprise.fit.coarse.comm <- functionFit(surprise.comm$surpriseIndex[surprise.comm$grain == "c" & surprise.comm$testCond == 'sample'], surprise.comm$sampleSize[surprise.comm$grain == "c" & surprise.comm$testCond == 'sample'], sampleSize)
surprise.fit.fine.comm <- functionFit(surprise.comm$surpriseIndex[surprise.comm$grain == "f" & surprise.comm$testCond == 'sample'], surprise.comm$sampleSize[surprise.comm$grain == "f" & surprise.comm$testCond == 'sample'], sampleSize)

surprise.functionFit <- surprise.fit.coarse.comm[[1]] %>% mutate(grain = "c")
surprise.functionFit <- rbind(surprise.functionFit, surprise.fit.fine.comm[[1]] %>% mutate(grain = "f"))
surprise.functionFit$grain <- as.factor(surprise.functionFit$grain)

#surprise.functionFit.comm.figure <- plotFunctionFits(surprise.functionFit, surprise.comm[surprise.comm$testCond == "sample",], "surpriseIndex", "Surprise Index (raw)")

for(i in 1:nrow(surprise.comm)){
  if(surprise.comm$grain[i] == 'c'){
    surprise.comm$elbow[i] = surprise.fit.coarse.comm[[3]][1]
  } else {
    surprise.comm$elbow[i] = surprise.fit.fine.comm[[3]][1]
  }
}

##everyday Activities 
surprise.fit.coarse.evAct <- functionFit(surprise.evAct$surpriseIndex[surprise.evAct$grain == "c" & surprise.evAct$testCond == "sample"], surprise.evAct$sampleSize[surprise.evAct$grain == 'c' & surprise.evAct$testCond == "sample"], sampleSize)
surprise.fit.fine.evAct <- functionFit(surprise.evAct$surpriseIndex[surprise.evAct$grain == "f" & surprise.evAct$testCond == "sample"], surprise.evAct$sampleSize[surprise.evAct$grain == 'f' & surprise.evAct$testCond == "sample"], sampleSize)

surprise.evAct.functionFit <- surprise.fit.coarse.evAct[[1]] %>% mutate(grain = "c")
surprise.evAct.functionFit <- rbind(surprise.evAct.functionFit, surprise.fit.fine.evAct[[1]] %>% mutate(grain = "f"))
surprise.evAct.functionFit$grain <- as.factor(surprise.evAct.functionFit$grain)

#surprise.functionFit.evAct.figure <- plotFunctionFits(surprise.evAct.functionFit, surprise.evAct[surprise.evAct$testCond == "sample",], "surpriseIndex", "Surprise Index (raw)") + ylim(0,15)

for(i in 1:nrow(surprise.evAct)){
  if(surprise.evAct$grain[i] == 'c'){
    surprise.evAct$elbow[i] = surprise.fit.coarse.evAct[[3]][1]
  } else {
    surprise.evAct$elbow[i] = surprise.fit.fine.evAct[[3]][1]
  }
}

Plot log10 Surprise Index over sample size

Log10 transformed peak-to-peak distance over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

surprise.comm.coarse.plot <- plotIndivAgreement(surprise.comm[surprise.comm$grain == 'c',], "logSurprise", "Surprise Index (log)", -1, 2, "Coarse")
surprise.comm.fine.plot <- plotIndivAgreement(surprise.comm[surprise.comm$grain == 'f',], "logSurprise", "Surprise Index (log)", -1, 2, "Fine")

surprise.evAct.coarse.plot <- plotIndivAgreement(surprise.evAct[surprise.evAct$grain == 'c', ],  "logSurprise", "Surprise Index (log)", -1, 2, "")
surprise.evAct.fine.plot <- plotIndivAgreement(surprise.evAct[surprise.evAct$grain == 'f', ],  "logSurprise", "Surprise Index (log)", -1, 2, "")

surprise.plot.combined <- ggarrange(surprise.comm.coarse.plot, surprise.comm.fine.plot,
                                      surprise.evAct.coarse.plot, surprise.evAct.fine.plot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

surprise.plot.combined

setwd('../plots')
ggsave("esMethod_surprise.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent", dpi = 300)

Plot Surprise Index function fits

Peak-to-peak values fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

surprise.comm.coarse.funcPlot <- plotFunctionFits(surprise.functionFit[surprise.functionFit$grain == 'c',], surprise.comm[surprise.comm$grain == 'c' & surprise.comm$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)
surprise.comm.fine.funcPlot <- plotFunctionFits(surprise.functionFit[surprise.functionFit$grain == 'f',], surprise.comm[surprise.comm$grain == 'f' & surprise.comm$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)

surprise.evAct.coarse.funcPlot <- plotFunctionFits(surprise.evAct.functionFit[surprise.evAct.functionFit$grain == 'c',], surprise.evAct[surprise.evAct$grain == 'c' & surprise.evAct$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)
surprise.evAct.fine.funcPlot <- plotFunctionFits(surprise.evAct.functionFit[surprise.evAct.functionFit$grain == 'f',], surprise.evAct[surprise.evAct$grain == 'f' & surprise.evAct$testCond == 'sample',], "surpriseIndex", "Surprise Index", 0, 15)

surprise.funcFit.plot.combined <- ggarrange(surprise.comm.coarse.funcPlot, surprise.comm.fine.funcPlot,
                                      surprise.evAct.coarse.funcPlot, surprise.evAct.fine.funcPlot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

surprise.funcFit.plot.combined

setwd('../plots')
ggsave("esMethod_surprise_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")

Surprise Index linear models and analysis

#Commercial 
surprise.comm.lmer <- lmer(surpriseIndex~sampleSize*grain*testCond + (1|movieName), data = surprise.comm)
#emm_options(pbkrtest.limit = 3000)
surprise.comm.poly <- summary(contrast(emmeans(surprise.comm.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
surprise.comm.rand_vs_samp <- summary(emmeans(surprise.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample 
surprise.comm.crossMov_vs_samp <- summary(emmeans(surprise.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample 
surprise.comm.crossMov_vs_rand <- summary(emmeans(surprise.comm.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 1, 0))), infer = T)# overall crossMov - random 

surprise.comm.pairwise <- summary(emmeans(surprise.comm.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

surprise.evAct.lmer <- lmer(surpriseIndex~sampleSize*grain*testCond + (1|movieName), data = surprise.evAct)
#emm_options(lmerTest.limit = 38400)
surprise.evAct.poly <- summary(contrast(emmeans(surprise.evAct.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
surprise.evAct.rand_vs_samp <- summary(emmeans(surprise.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample 
surprise.evAct.crossMov_vs_samp <- summary(emmeans(surprise.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample 
surprise.evAct.crossMov_vs_rand <- summary(emmeans(surprise.evAct.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random

surprise.evAct.pairwise <- summary(emmeans(surprise.evAct.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

Surprise Index random - sample contrast commercial-lab

surprise.comm.rand_vs_samp$contrasts
## grain = c:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       3.434135 0.03548559 Inf  3.364585  3.503686  96.775  <.0001
## 
## grain = f:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       6.705420 0.03548559 Inf  6.635870  6.774971 188.962  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Surprise Index random - sample contrast everyday-online

surprise.evAct.rand_vs_samp$contrasts
## grain = c:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       2.549908 0.02249992 Inf  2.505809  2.594007 113.330  <.0001
## 
## grain = f:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       4.052439 0.02249992 Inf  4.008340  4.096538 180.109  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Surprise Index cross movie - sample contrast commercial-lab

surprise.comm.crossMov_vs_samp$contrasts
## grain = c:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       3.570111 0.03548559 Inf  3.500561  3.639662 100.607  <.0001
## 
## grain = f:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       6.679965 0.03548559 Inf  6.610415  6.749516 188.244  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Surprise Index cross movie - sample contrast everyday-online

surprise.evAct.crossMov_vs_samp$contrasts
## grain = c:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       2.262476 0.02249992 Inf  2.218377  2.306575 100.555  <.0001
## 
## grain = f:
##  contrast estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       3.885204 0.02249992 Inf  3.841105  3.929303 172.676  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Surprise Index cross movie - random contrast commercial-lab

surprise.comm.crossMov_vs_rand$contrasts
## grain = c:
##  contrast    estimate         SE  df   asymp.LCL  asymp.UCL z.ratio
##  c1        0.13597635 0.03548559 Inf  0.06642588 0.20552683   3.832
##  p.value
##   0.0001
## 
## grain = f:
##  contrast    estimate         SE  df   asymp.LCL  asymp.UCL z.ratio
##  c1       -0.02545493 0.03548559 Inf -0.09500541 0.04409554  -0.717
##  p.value
##   0.4732
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Surprise Index cross movie - random contrast everyday-online

surprise.evAct.crossMov_vs_rand$contrasts
## grain = c:
##  contrast  estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.2874323 0.02249992 Inf 0.2433333 0.3315313  12.775  <.0001
## 
## grain = f:
##  contrast  estimate         SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.1672345 0.02249992 Inf 0.1231355 0.2113335   7.433  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Surprise Index smallest sample size for significant pairwise effect commercial-lab

surprise.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 4 x 10
## # Groups:   grain, contrast [4]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c         1.77 0.142   Inf      1.49      2.04
## 2 sample … 2          c         1.60 0.142   Inf      1.33      1.88
## 3 sample … 2          f         2.33 0.142   Inf      2.05      2.60
## 4 sample … 2          f         2.32 0.142   Inf      2.04      2.60
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Surprise Index smallest sample size for significant pairwise effect everyday-online

surprise.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 4 x 10
## # Groups:   grain, contrast [4]
##   contrast sampleSize grain estimate     SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>  <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c         1.28 0.0900   Inf      1.11      1.46
## 2 sample … 2          c         1.21 0.0900   Inf      1.03      1.38
## 3 sample … 2          f         1.88 0.0900   Inf      1.70      2.05
## 4 sample … 2          f         1.91 0.0900   Inf      1.73      2.09
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Surprise Index smallest pairwise effect commercial-lab

surprise.comm.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 4 x 10
## # Groups:   grain, contrast [4]
##   contrast sampleSize grain estimate    SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl> <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 4          c         1.44 0.142   Inf      1.16      1.72
## 2 sample … 4          c         1.33 0.142   Inf      1.05      1.61
## 3 sample … 2          f         2.33 0.142   Inf      2.05      2.60
## 4 sample … 2          f         2.32 0.142   Inf      2.04      2.60
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Surprise Index smallest pairwise effect everyday-online

surprise.evAct.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 4 x 10
## # Groups:   grain, contrast [4]
##   contrast sampleSize grain estimate     SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>  <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c         1.28 0.0900   Inf      1.11      1.46
## 2 sample … 2          c         1.21 0.0900   Inf      1.03      1.38
## 3 sample … 2          f         1.88 0.0900   Inf      1.70      2.05
## 4 sample … 2          f         1.91 0.0900   Inf      1.73      2.09
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Surprise Index polynomial fit commercial-lab

surprise.comm.poly[surprise.comm.poly$testCond == "sample",]
##     contrast grain testCond   estimate         SE  df asymp.LCL  asymp.UCL
## 13    linear     c   sample  157.21190   3.701404 Inf  149.9573  164.46652
## 14 quadratic     c   sample  -43.69346   7.585620 Inf  -58.5610  -28.82591
## 15     cubic     c   sample -392.21574 100.757076 Inf -589.6960 -194.73550
## 16    linear     f   sample  287.41961   3.701404 Inf  280.1650  294.67423
## 17 quadratic     f   sample -211.13912   7.585620 Inf -226.0067 -196.27157
## 18     cubic     f   sample  351.69993 100.757076 Inf  154.2197  549.18017
##       z.ratio       p.value
## 13  42.473591  0.000000e+00
## 14  -5.760038  8.409524e-09
## 15  -3.892687  9.914009e-05
## 16  77.651517  0.000000e+00
## 17 -27.834127 1.676671e-170
## 18   3.490573  4.819859e-04

Surprise Index polynomial fit everyday-online

surprise.evAct.poly[surprise.evAct.poly$testCond == "sample",]
##     contrast grain testCond  estimate        SE  df  asymp.LCL asymp.UCL
## 13    linear     c   sample  78.15937  2.346904 Inf   73.55952  82.75922
## 14 quadratic     c   sample -31.27365  4.809722 Inf  -40.70053 -21.84676
## 15     cubic     c   sample -74.87400 63.885815 Inf -200.08790  50.33990
## 16    linear     f   sample  97.38823  2.346904 Inf   92.78838 101.98807
## 17 quadratic     f   sample -85.25332  4.809722 Inf  -94.68021 -75.82644
## 18     cubic     f   sample 577.77115 63.885815 Inf  452.55725 702.98505
##       z.ratio       p.value
## 13  33.303180 3.472692e-243
## 14  -6.502173  7.916788e-11
## 15  -1.171997  2.411982e-01
## 16  41.496467  0.000000e+00
## 17 -17.725208  2.679205e-70
## 18   9.043810  1.513047e-19

Predictive value and detection accuracy

Load dataset

setwd('../data/bootstrapped')
normativeHit.comm <- read.delim("esMethods_normativeHitRate_bin_100Iterations_commercial.txt", head = TRUE)

predictive.value<- normativeHit.comm %>% dplyr::select(Iteration, sampleSize, movieName, grain, testCond, normativePPV, normativeFalsePPV)
predictive.value$dPrime <- qnorm(predictive.value$normativePPV) - qnorm(predictive.value$normativeFalsePPV)
predictive.value <- predictive.value %>% dplyr::filter(is.finite(predictive.value$dPrime))
predictive.value$sampleSize <- as.factor(predictive.value$sampleSize)

detection.sensitivity <- normativeHit.comm %>% dplyr::select(Iteration, sampleSize, movieName, grain, testCond, normativeHit, normativeFA, normativeMiss, normativeCR)
detection.sensitivity$dPrime <- qnorm(detection.sensitivity$normativeHit) - qnorm(detection.sensitivity$normativeFA)
detection.sensitivity <- detection.sensitivity %>% dplyr::filter(is.finite(detection.sensitivity$dPrime))
detection.sensitivity$sampleSize <- as.factor(detection.sensitivity$sampleSize)

Fit predictive value and detection accuracy values to growth and decay functions

#zalladPrime
predictive.fit.coarse <- functionFit(predictive.value$dPrime[predictive.value$grain == "c" & predictive.value$testCond == 'sample'], predictive.value$sampleSize[predictive.value$grain == "c" & predictive.value$testCond == 'sample'], sampleSize)
## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
##   non-finite value supplied by optim
## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
##   non-finite value supplied by optim
predictive.fit.coarse <- functionFit(predictive.value$dPrime[predictive.value$grain == "f" & predictive.value$testCond == 'sample'], predictive.value$sampleSize[predictive.value$grain == "f" & predictive.value$testCond == 'sample'], sampleSize)
## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
##   non-finite value supplied by optim
## Error in optim(startVec, opfct, hessian = TRUE, method = optMethod, control = list(maxit = maxIt,  : 
##   non-finite value supplied by optim
predictive.functionFit <- predictive.fit.coarse[[1]] %>% mutate(grain = "c")
predictive.functionFit <- rbind(predictive.functionFit, predictive.fit.coarse[[1]] %>% mutate(grain = "f"))
predictive.functionFit$grain <- as.factor(predictive.functionFit$grain)

for(i in 1:nrow(predictive.value)){
  if(predictive.value$grain[i] == 'c'){
    predictive.value$elbow[i] = predictive.fit.coarse[[3]][1]
  } else {
    predictive.value$elbow[i] = predictive.fit.coarse[[3]][1]
  }
}


##commercial 
detection.sensitivity.fit.coarse.comm <- functionFit(detection.sensitivity$dPrime[detection.sensitivity$grain == "c" & detection.sensitivity$testCond == 'sample'], detection.sensitivity$sampleSize[detection.sensitivity$grain == "c" & detection.sensitivity$testCond == 'sample'], sampleSize)
detection.sensitivity.fit.fine.comm <- functionFit(detection.sensitivity$dPrime[detection.sensitivity$grain == "f" & detection.sensitivity$testCond == 'sample'], detection.sensitivity$sampleSize[detection.sensitivity$grain == "f" & detection.sensitivity$testCond == 'sample'], sampleSize)

detection.sensitivity.functionFit <- detection.sensitivity.fit.coarse.comm[[1]] %>% mutate(grain = "c")
detection.sensitivity.functionFit <- rbind(detection.sensitivity.functionFit, detection.sensitivity.fit.fine.comm[[1]] %>% mutate(grain = "f"))
detection.sensitivity.functionFit$grain <- as.factor(detection.sensitivity.functionFit$grain)


for(i in 1:nrow(detection.sensitivity)){
  if(detection.sensitivity$grain[i] == 'c'){
    detection.sensitivity$elbow[i] = detection.sensitivity.fit.coarse.comm[[3]][1]
  } else {
    detection.sensitivity$elbow[i] = detection.sensitivity.fit.fine.comm[[3]][1]
  }
}

Plot predictive value and detection accuracy over sample size

Predictive value and detection accuracy over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

predictive.value.comm.coarse.plot <- plotIndivAgreement(predictive.value[predictive.value$grain == 'c',], "dPrime", "Predictive value", -6, 1, "Coarse")
predictive.value.comm.fine.plot <- plotIndivAgreement(predictive.value[predictive.value$grain == 'f',],"dPrime", "Predictive value", -6, 1, "Fine")

predictive.value.plot <- ggarrange(predictive.value.comm.coarse.plot, predictive.value.comm.fine.plot,
                                       nrow = 1, ncol = 2,
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

detection.sensitivity.coarse.plot <- plotIndivAgreement(detection.sensitivity[detection.sensitivity$grain == 'c',], "dPrime", "Detection accuracy", -1, 2, "Coarse")
detection.sensitivity.fine.plot <- plotIndivAgreement(detection.sensitivity[detection.sensitivity$grain == 'f',], "dPrime", "Detection accuracy", -1, 2, "Fine")

detection.sensitivity.plot <- ggarrange(detection.sensitivity.coarse.plot, detection.sensitivity.fine.plot,
                                       nrow = 1, ncol = 2,
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

predictive.detection.accuracy <- ggarrange(predictive.value.plot, detection.sensitivity.plot,
                                  nrow = 2, ncol = 1,
                                  common.legend = TRUE, legend = "bottom",
                                  labels = c("A","B"),
                                  align = "hv")

predictive.detection.accuracy

setwd('../plots')
ggsave("esMethod_predictive&detection.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")

Plot predictive value and detection accuracy function fits

Predictive value and detection accuracy fit to various growth and decay functions over increasing sample sizes for segmentation of: A) commercial-lab and B) everyday-online.

predictive.coarse.funcPlot <- plotFunctionFits(predictive.functionFit[predictive.functionFit$grain == 'c',], predictive.value[predictive.value$grain == 'c' & predictive.value$testCond == 'sample',], "dPrime", "Predictive value", -4, 0)
predictive.fine.funcPlot <- plotFunctionFits(predictive.functionFit[predictive.functionFit$grain == 'f',], predictive.value[predictive.value$grain == 'f' & predictive.value$testCond == 'sample',], "dPrime", "Predictive value", -4, 0)

detection.coarse.funcPlot <- plotFunctionFits(detection.sensitivity.functionFit[detection.sensitivity.functionFit$grain == 'c',], detection.sensitivity[detection.sensitivity$grain == 'c' & detection.sensitivity$testCond == 'sample',], "dPrime", "Detection accuracy", -1, 2)
detection.fine.funcPlot <- plotFunctionFits(detection.sensitivity.functionFit[detection.sensitivity.functionFit$grain == 'f',], detection.sensitivity[detection.sensitivity$grain == 'f' & detection.sensitivity$testCond == 'sample',], "dPrime", "Detection accuracy", -1, 2)


predictive.detection.funcFit.plot.combined <- ggarrange(predictive.coarse.funcPlot, predictive.fine.funcPlot,
                                      detection.coarse.funcPlot, detection.fine.funcPlot,
                                       nrow = 2, ncol = 2,
                                       labels = c("A","", "B", ""),
                                       common.legend = TRUE, legend = "bottom",
                                       align = "hv")

predictive.detection.funcFit.plot.combined

setwd('../plots')
ggsave("esMethod_predictive&detection_functionFit.pdf", device = cairo_pdf, width = 25, height = 19, units = "cm", bg = "transparent")

Predictive value and detection accuracy linear models and analysis

#Predictive value  
predictive.value.lmer <- lmer(dPrime~sampleSize*grain*testCond + (1|movieName), data = predictive.value)
#emm_options(pbkrtest.limit = 3000)

predictive.rand_vs_samp <- summary(emmeans(predictive.value.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample 
predictive.crossMov_vs_samp <- summary(emmeans(predictive.value.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample 
predictive.crossMov_vs_rand <- summary(emmeans(predictive.value.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random 

predictive.poly <- summary(contrast(emmeans(predictive.value.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast

predictive.pairwise <- summary(emmeans(predictive.value.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

#Detection sensitivity  
detection.sensitivity.lmer <- lmer(dPrime~sampleSize*grain*testCond + (1|movieName), data = detection.sensitivity)
#emm_options(pbkrtest.limit = 19029)
detection.poly <- summary(contrast(emmeans(detection.sensitivity.lmer,'sampleSize', by = c( 'grain', 'testCond')), 'poly', max.degree = 3), infer = TRUE)#Polynomial contrast
detection.rand_vs_samp <- summary(emmeans(detection.sensitivity.lmer, c("testCond"), by = "grain", contr = list(c1 = c(0, -1, 1))), infer = T)# overall random - sample 
detection.crossMov_vs_samp <- summary(emmeans(detection.sensitivity.lmer, c("testCond"), by = "grain", contr = list(c1 = c(-1, 0, 1))), infer = T)# overall crossMov - sample 
detection.crossMov_vs_rand <- summary(emmeans(detection.sensitivity.lmer, c("testCond"), by = "grain", contr = list(c1 = c(1, -1, 0))), infer = T)# overall crossMov - random 

detection.pairwise <- summary(emmeans(detection.sensitivity.lmer, "testCond",by = c("sampleSize", "grain"), weights = "proportional", contr = "revpairwise", adjust = "none")$contrasts, infer = TRUE)#random vs sample pairwise/ sample size

Predictive value random - sample contrast commercial-lab

predictive.rand_vs_samp$contrasts
## grain = c:
##  contrast estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       1.585063 0.009201382 Inf  1.567029  1.603097 172.264  <.0001
## 
## grain = f:
##  contrast estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       1.115986 0.009081682 Inf  1.098186  1.133786 122.883  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Detection value random - sample contrast commercial-lab

detection.rand_vs_samp$contrasts
## grain = c:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.8096764 0.004325468 Inf 0.8011986 0.8181541 187.188  <.0001
## 
## grain = f:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.6529486 0.004269173 Inf 0.6445812 0.6613161 152.945  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Predictive value cross movie - sample contrast commercial-lab

predictive.crossMov_vs_samp$contrasts
## grain = c:
##  contrast estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       1.481705 0.009156284 Inf  1.463759  1.499651 161.824  <.0001
## 
## grain = f:
##  contrast estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       1.125738 0.009081682 Inf  1.107938  1.143538 123.957  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Detection value cross movie - sample contrast commercial-lab

detection.crossMov_vs_samp$contrasts
## grain = c:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.7738612 0.004304249 Inf 0.7654250 0.7822974 179.790  <.0001
## 
## grain = f:
##  contrast  estimate          SE  df asymp.LCL asymp.UCL z.ratio p.value
##  c1       0.6622736 0.004269173 Inf 0.6539062 0.6706410 155.129  <.0001
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Predictive value cross movie - random contrast commercial-lab

predictive.crossMov_vs_rand$contrasts
## grain = c:
##  contrast     estimate          SE  df   asymp.LCL   asymp.UCL z.ratio
##  c1        0.103357684 0.009176672 Inf  0.08537174 0.121343630  11.263
##  p.value
##   <.0001
## 
## grain = f:
##  contrast     estimate          SE  df   asymp.LCL   asymp.UCL z.ratio
##  c1       -0.009752199 0.009079522 Inf -0.02754774 0.008043337  -1.074
##  p.value
##   0.2828
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Detection value cross movie - random contrast commercial-lab

detection.crossMov_vs_rand$contrasts
## grain = c:
##  contrast    estimate          SE  df   asymp.LCL     asymp.UCL z.ratio
##  c1        0.03581518 0.004313832 Inf  0.02736022  0.0442701338   8.302
##  p.value
##   <.0001
## 
## grain = f:
##  contrast    estimate          SE  df   asymp.LCL     asymp.UCL z.ratio
##  c1       -0.00932495 0.004268157 Inf -0.01769038 -0.0009595155  -2.185
##  p.value
##   0.0289
## 
## Results are averaged over the levels of: sampleSize 
## Confidence level used: 0.95

Predictive value smallest sample size for significant pairwise effect commercial-lab

predictive.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 5 x 10
## # Groups:   grain, contrast [5]
##   contrast sampleSize grain estimate     SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>  <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.737  0.0398   Inf    0.659      0.815
## 2 sample … 2          c       0.914  0.0422   Inf    0.831      0.997
## 3 random … 4          c       0.0978 0.0377   Inf    0.0240     0.172
## 4 sample … 2          f       0.439  0.0365   Inf    0.368      0.511
## 5 sample … 2          f       0.516  0.0365   Inf    0.444      0.587
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Detection value smallest sample size for significant pairwise effect commercial-lab

detection.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(sampleSize == min(as.numeric(as.character(sampleSize))))
## # A tibble: 6 x 10
## # Groups:   grain, contrast [6]
##   contrast sampleSize grain estimate     SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>  <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.343  0.0187   Inf  0.306       0.380 
## 2 sample … 2          c       0.501  0.0199   Inf  0.462       0.540 
## 3 random … 4          c       0.0399 0.0177   Inf  0.00523     0.0746
## 4 sample … 2          f       0.254  0.0171   Inf  0.220       0.287 
## 5 sample … 2          f       0.311  0.0171   Inf  0.278       0.345 
## 6 random … 10         f       0.0340 0.0171   Inf  0.000543    0.0675
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Predictive value smallest pairwise effect commercial-lab

predictive.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 5 x 10
## # Groups:   grain, contrast [5]
##   contrast sampleSize grain estimate     SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>  <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.737  0.0398   Inf    0.659      0.815
## 2 sample … 2          c       0.914  0.0422   Inf    0.831      0.997
## 3 random … 4          c       0.0978 0.0377   Inf    0.0240     0.172
## 4 sample … 2          f       0.439  0.0365   Inf    0.368      0.511
## 5 sample … 2          f       0.516  0.0365   Inf    0.444      0.587
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Detection value smallest pairwise effect commercial-lab

detection.pairwise %>% dplyr::filter(p.value < .05, z.ratio > 0) %>% dplyr::group_by(grain, contrast) %>% filter(z.ratio == min(z.ratio))
## # A tibble: 6 x 10
## # Groups:   grain, contrast [6]
##   contrast sampleSize grain estimate     SE    df asymp.LCL asymp.UCL
##   <fct>    <fct>      <fct>    <dbl>  <dbl> <dbl>     <dbl>     <dbl>
## 1 sample … 2          c       0.343  0.0187   Inf  0.306       0.380 
## 2 sample … 2          c       0.501  0.0199   Inf  0.462       0.540 
## 3 random … 4          c       0.0399 0.0177   Inf  0.00523     0.0746
## 4 sample … 2          f       0.254  0.0171   Inf  0.220       0.287 
## 5 sample … 2          f       0.311  0.0171   Inf  0.278       0.345 
## 6 random … 10         f       0.0340 0.0171   Inf  0.000543    0.0675
## # … with 2 more variables: z.ratio <dbl>, p.value <dbl>

Predictive value polynomial fit commercial-lab

predictive.poly[predictive.poly$testCond == "sample",]
##     contrast grain testCond  estimate         SE  df asymp.LCL asymp.UCL
## 13    linear     c   sample  33.58800  0.9738438 Inf  31.67930  35.49670
## 14 quadratic     c   sample -26.37411  2.0076247 Inf -30.30898 -22.43924
## 15     cubic     c   sample  90.12431 26.5866237 Inf  38.01549 142.23314
## 16    linear     f   sample  27.82441  0.9482521 Inf  25.96587  29.68295
## 17 quadratic     f   sample -30.80938  1.9440620 Inf -34.61968 -26.99909
## 18     cubic     f   sample 237.24458 25.8205091 Inf 186.63731 287.85184
##       z.ratio       p.value
## 13  34.490135 1.127637e-260
## 14 -13.136972  2.021558e-39
## 15   3.389837  6.993427e-04
## 16  29.342840 2.948948e-189
## 17 -15.847943  1.452581e-56
## 18   9.188222  3.993872e-20

Detection value polynomial fit commercial-lab

detection.poly[detection.poly$testCond == "sample",]
##     contrast grain testCond  estimate         SE  df asymp.LCL asymp.UCL
## 13    linear     c   sample  16.80430  0.4577924 Inf  15.90704  17.70156
## 14 quadratic     c   sample -13.06385  0.9437606 Inf -14.91359 -11.21411
## 15     cubic     c   sample  43.90222 12.4980363 Inf  19.40652  68.39792
## 16    linear     f   sample  15.92080  0.4457602 Inf  15.04713  16.79448
## 17 quadratic     f   sample -17.42398  0.9138766 Inf -19.21514 -15.63281
## 18     cubic     f   sample 127.27518 12.1378629 Inf 103.48541 151.06495
##       z.ratio       p.value
## 13  36.707249 5.595930e-295
## 14 -13.842334  1.415362e-43
## 15   3.512729  4.435288e-04
## 16  35.716072 2.225900e-279
## 17 -19.066008  4.838610e-81
## 18  10.485798  1.003944e-25

#Extra analyses

setwd('../data/raw')
ts.data.commercial <- read.delim("esMethods_Commercial_TimeSeries.txt", head = TRUE)
ts.data.evAct <- read.delim("esMethods_EvAct_TimeSeries.txt", head = TRUE)
names(ts.data.commercial)[names(ts.data.commercial) == "movName"] <- "movieName"
names(ts.data.evAct)[names(ts.data.evAct) == "movName"] <- "movieName"

Correlate correlation between segmentation of different movies

comm.gp.bp <- ts.data.commercial %>% dplyr::group_by(grain, movieName, timeSeries) %>% dplyr::summarise(bp = mean(subject.bp))
evAct.gp.bp <- ts.data.evAct %>% dplyr::group_by(grain, movieName, timeSeries) %>% dplyr::summarise(bp = mean(subject.bp))

comm.corr.fine <- cor.test(comm.gp.bp$bp[comm.gp.bp$movieName == "3Iron" & comm.gp.bp$grain == "f"],
                           comm.gp.bp$bp[comm.gp.bp$movieName == "Corn" & comm.gp.bp$grain == "f"])

comm.corr.coarse <- cor.test(comm.gp.bp$bp[comm.gp.bp$movieName == "3Iron" & comm.gp.bp$grain == "c"],
                           comm.gp.bp$bp[comm.gp.bp$movieName == "Corn" & comm.gp.bp$grain == "c"])

evAct.corr.fine.bed_coffee <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "f"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.fine.bed_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "f"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.fine.bed_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "f"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "f"])

evAct.corr.coarse.bed_coffee <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "c"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.coarse.bed_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "c"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.coarse.bed_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_bed.mp4" & evAct.gp.bp$grain == "c"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "c"])


evAct.corr.fine.coffee_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "f"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "f"])
evAct.corr.fine.coffee_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "f"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "f"])


evAct.corr.coarse.coffee_dishes <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "c"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "c"])
evAct.corr.coarse.coffee_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_coffee.mp4" & evAct.gp.bp$grain == "c"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "c"])


evAct.corr.fine.dishes_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "f"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "f"])

evAct.corr.coarse.dishes_laundry <- cor.test(evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_dishes.mp4" & evAct.gp.bp$grain == "c"],
                           evAct.gp.bp$bp[evAct.gp.bp$movieName == "w_laundry.mp4" & evAct.gp.bp$grain == "c"])